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A detailed modelling of the chemically rich clumps along 
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ABSTRACT 

In order to investigate the origin and the structure of the low velocity, chemically rich 
clumps observed along the lobes of low- and intermediate-mass outflows, we construct 
a detailed model of the SI clump along the CBS outflow. We use a time-dependent 
chemical model coupled with a radiative transfer model to reproduce the observed 
line profile for a direct comparison with previous observations of this clump. We find 
that the simultaneous fitting of multiple species and transitions is a powerful tool 
in constraining the physical parameters of the gas. Different scenarios for the clump 
formation have been investigated. The models that better reproduce all the observed 
lines are those where the clump is formed, at least partially, before the advent of the 
outflow; with the advent of the outflow the clump undergoes a short period of non- 
dissociative shock and the consequent release of the icy mantle together with the high 
temperature chemistry leads to the observed chemical enrichment. Our results also 
suggest the presence of substructure within the clump: a more extended component 
traced by CS, SO and the lower energy transitions (3k-2k and 2k-1k) of CH3OH, and 
a more compact component traced by H2CO, SO2 and the higher energy transitions 
(5k-4k) of CH3OH. 

Key words: ISM: individual: CB3 - ISM: jets and outflows - ISM: molecules. 



1 INTRODUCTION 

Bipolar molecular outflows are ubiquitously present around 
Young Stellar Objects (YSOs) and are associated with th e 
very early stages of star formation (e.g. lRicher et alj|200(j) . 
A large range of excitation conditions and a significant 
chemical differentiation are present along the lobes of the 
outflows. The shock-excited gas cools mainly radiatively via 
the emission of atomic and molecular species whose abun- 
dances are strongly enhanced by the shocks generated by the 
interaction between the protostellar wind and the molec- 
ular cloud. High velocity (v^lOO km s~^) clumps, the so 
called molecular bullets, are detected along the axis of some 
outflows. The prototype of such bullets are those observed 
in Si O in the lobes of the L1148 outflow ( jBachiller et al . 
llQQlT) but they generally show weak emission of other molec- 
ular lines. A few outflows however, e.g. LI 15 7 and CBS 
jBachiUer et alj200ll . lCodefla fc Bachilleill99gll . have a par- 
ticular rich emission spectrum with molecular clumps at rel- 
atively low velocity (v5C 10 km s^^). These molecular clumps 
show emission from molecular species usually unobserved 
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along outflows and they have dimensions which do not ex- 
ceed 0.1 pc. While the high velocities and high abundance 
of SiO suggest that the molecular bullets are most likely as- 
sociated wit h mini-bow shocks fo rmed by the outflow prop- 
agation fe.g. iDeutrev et al]|l997^ . the origin of t he low ve- 
locity clumps is not yet clear. In a recent study IViti et alJ 
(J2004) (thereafter Paper I), have investigated two possibili- 
ties for the origin of the low velocity clumps: i) the clumps 
are formed by episodic mass loss of the forming object or 
ii) the clumps are density structure pre-existing the advent 
of the outflow and the effect of the outflow is to acceler- 
ate and shock them, altering their chemistry. Paper I used a 
chemical model to simulate the clump formation and its sub- 
sequent interaction with the outflow and inferred that the 
low velocity outflow clumps are probably density structures 
formed, at least partly, prior of the advent of the outflow. 
With the advent of the outflow, not only the temperature of 
the clumps increases (~ 100 K), but the clumps also undergo 
a short period of non-dissociative shock. This preliminary 
conclusion was reached by comparing the column densities 
derived by the chemical models with existing single dis h ob- 
servations alon g the CB3 iCodell a fc BachilleJ Il999l) and 
L1157 outflows JBachiller et alJl2001f) . The shortcomings of 
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such an analysis were that the estimations of the column 
densities derived from the chemical models, as well as those 
derived from the observations, suffer from high uncertain- 
ties due to the fact that some of the parameters required 
to calculate them may be unknown. In particular, to esti- 
mate the observed column densities arbitrary assumptions 
have to be made regarding LTE conditions, excitation tem- 
peratures and the lines being optically thin. On the other 
hand, the theoretical calculation of column densities from 
the chemical model required knowledge of the geometry of 
the emitting region. However, when dealing with small emit- 
ting regions, such as the low velocity clumps along molecular 
outflows, the observations (usually single dish with spatial 
resolution of few tens of arcsec) do not constrain such ge- 
ometry. If more transitions of the same species are available, 
one can use a Large Velocity Gradient (LVG) model to de- 
rive the density and temperature of the emitting gas but 
even in such cases often the results have a large uncertainty 
because the observed line ratios are not always sensitive to 
the temperature or density. 

One way to reduce the number of arbitrary assumptions 
is to model directly the observed line emission using a radia- 
tive transfer model that predicts the l ine profile, a method 
already adopted by o ther authors (e.g. |j 0rgesen et alj|2004l 
iDotv et al.ll2004l and ISchoier et aljl2002l) . In this work we 
attempt such an approach by coupling the time dependent 
chemical model of Paper I with the radiative transfer mo del 
Spherical Multi-Mol (SMMOL) jRawlings fc Yatej200lH . In 
particular, we expand the grids of chemical models com- 
puted in Paper I and use them as input parameters for the 
radiative transfer model to simulate the line profiles of sev- 
eral species and, where possible, of multiple transitions of 
the same species. We then simultaneously fit such profiles 
with the observations of CBS and attempt to constrain the 
physical parameters of the emitting gas and thus their for- 
mation scenario. The paper is organized as follow: in Section 
2 both the chemical and the radiative transfer models are 
briefly described. In Section 3 the models are compared with 
the lines detected towards a clump in the CBS outflow. The 
results are presented in Section 4 and discussed in Section 
5. 



2 THE CHEMICAL AND RADIATIVE 
TRANSFER MODELS 

A flow diagram of the method used to compare the ob- 
served data with the theoretical models is shown in Fig. 
A first order estimate of the general physical conditions 
of the clump, such as the size, the density and the temper- 
ature, is derived from the observations and used to describe 
the physics of the phenomenon in the chemical model. The 
chemical model predicts the fractional abundances; these 
are given as input to the radiative transfer model together 
with the density, the temperature and the beam of the tele- 
scope to produce the line profiles that are directly compared 
with the observed profiles. Recently, the same method has 
been used for modelling the envelope of low mass protostars 
Jjorgesen et al J 120041: 120051: ISchoier et a'i]l2002l : iDotv et alJ 
|200j) and it was found to be a powerful tool to probe the 
physical and chemical structure of the protostellar envelopes 
as well as the age of the protostar. 
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Figure 1. A flow diagram of the scheme applied to compare the 
observed data to the theoretical models in order to derive the best 
fit model. 



The time-dependent chemical model is described in de- 
tail in Paper I. Here we summarize the important charac- 
teristics: the chemical network is taken from the UMIST 
database and includes 221 species involved in 3194 gas-phase 
and grain reactions. The model is a two-phase calculation. 
Phase I simulates the formation of a pre-existing clump, or 
of the homogeneous dark cloud, depending on what scenario 
we are modelling, and includes freeze-out of the gas on to 
the grains (see Paper I and below). Phase II simulates the 
advent of the outflow by warming the clump to 100-200 K 
or, in the presence of a non-dissociative shock, by increasing 
the gas temperature to 1000 K for 100 yr and then cool- 
ing it down to 100-200 K. As a reminder, the four scenarios 
considered here (as in Paper I) are: 

• Grid A: a pre-existent clump at uniform density warmed 
by the outflow. In this paper we introduce new models, called 
As, where the clump at uniform single density is shocked by 
the outflow; 

• Grid B: a pre-existing clump with a density structure 
shocked by the outflow; 

• Grid C: a pre-existing clump with a density structure 
warmed by the outflow; 

• Grid D: a clump formed and warmed or shocked by the 
outflow. 

For each grid difl'erent conditions are investigated, 
changing the model parameters: the depletion efficiency, the 
type of collapse (free fall or retarded) , the grain mantle evap- 
oration (instantaneous or time-dependent), the final density 
and temperature, the initial sulphur abundance (see Section 
4.1), the fraction of frozen H2CO and CO converted into 
CH3OH (see Section 4.4), for a total of SS models. Each 
model is time-dependent and evolves for 10^ yr. See Paper 
I for a full description of the model and its assumptions, as 
well as of each scenario. 

The chemical abundances produced by the chemical 
model together with the density and temperature profile 
are used as input for the radiative transfer model SMMOL 
JRawlings fc YateJ200ir) . an approximate A-iteration (ALI) 
code that solves multi-level non-LTE radiative transfer prob- 
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Table 1. List of the molecules and transitions observed in the 
SI clump of the CBS outflow that are modeled with the radiative 
transfer code. 



specie 



line 



u (GHz) 



CO 


2-1 


230.538 


cs 


3-2 


146.967 


so 


65 — 54 


219.949 


so 


43-32 


138.178 


S02 


3l3-202 


104.029 


H2C0 


321-220 


218.760 


CH30H 


5o-4o A+ 


241.791 


CH30H 


5_i-4_i E 


241.767 


CH30H 


5o-4o E 


241.700 


CH30H 


3o-2o A+ 


145.103 


CH30H 


3-i-2_i E 


145.097 


CH30H 


3o-2o E 


145.094 


CH30H 


20-I0A+ 


96.741 


CH30H 


2-1-1-1 E 


96.739 



lems. The molecular data needed for the c alculation are 
taken from the LAMBDA^ public database l|Schoier et alJ 
120051) . except the CH3OH-A data which we obtained from 
D. Flower (private communication). As first step the code 
calculates the total radiation field and the level populations 
assuming L TE and the interstellar radiation field as input 
continuum (|Blacyll99J). It then re-calculates the total ra- 
diation field, checks for convergence and repeats the process 
until convergence is achieved. At each radial point the code 
generates the level populations and the radiation field. The 
emergent intensity distributions are then convolved with the 
telescope beam, so that the model directly predicts the line 
profiles for a given source as observed with a given telescope 
(e.g. the IRAM 30-m telescope in this work). The radia- 
tive transfer model has been successfully bonchmarkcd with 
other radiative transfer codes (ivan Zadclhoff ot al. 2002). 
For each model the radiative transfer program has been run 
for a selected subset of the time steps of the chemical model 
for a total of 19 time steps starting from 10^ yr after the 
formation of the clump to 10^ yr. 



the radiative transfer code has been run. Most of these lines 
have a Gaussian profile with similar width. We adopt a mi- 
croturbulent velocity of 5 km s~^ for our radiative transfer 
modelling. This is equivalent to modelling the observed line 
profiles with a single gaussian line with a fixed FWHM of 10 
km s~^. We assume a clump size of 0.12 pc as derived from 
the IRAM maps. However, we note that most of the emis- 
sion from the clump is unresolved by the large IRAM-30m 
beam (HPBW = 10"- 29"), so 0.12 pc is most likely an over- 
estimate of the real size of the clump. In Table |5| we list the 
chemical models computed and their parameters: the name 
of the model (column 1); the presence (Y) or absence (N) of 
non-dissociative shock (column 2); the density of the clump 
(column 3), for models with a density profile we give the two 
extremes of the range; the final temperature (column 4); the 
percentage of gas depleted into grains at the end of the Phase 
I (column 5); the initial sulfur abundance (column 6); and 
any other relevant parameters that were varied (column 7) 
such as the type of collapse simulated in Phase I (free-fall 
or retarded), the percentage of CO or H2CO converted into 
CH3OH on the grain surfaces or the type of evaporation from 
the grains once the outflow ap proaches the clumps (ins tanta- 
neous or time-dependent as in lViti fc Williamsl il999l) ). For 
the models of Grid D we also give the velocity of the outflow. 
For this work, we extend the grids of models already pub- 
lished in Paper I. In particular, we computed new models 
for: clumps at uniform density that are shocked (Grid As); 
clumps with a density proflle from 5x10^ cm~^ to 10® cm~'^ 
at a temperature of T=200 K; clumps with a density profile 



cm 

-3 



to 5x10^ cm""* and from 3x10'' cm"'' 
From Paper I, we found that most of the 



from 2.5x10'' 

to 6x10^ cm' 

models from Grid D, where the clumps are formed by the 

compression of the incoming outfiow, give unphysical results 

and thus in this paper we restrict our analysis to models Dl, 

D3, D4, D4-shock and D8 only. 

To fit the lines listed in Table 1 we used a x^ method 
that simultaneously fits the total integrated fiux of the ob- 
served transitions of various species. The formula used is: 



X 



-y 



FluXmod{i) - FluXobs{i) 



FluXobsii) 



(1) 



3 MODELLING THE CLUMPS ALONG CBS 

CB3 is a Bok globule at a distance of 2500 pc associated with 
star formation. In particular, it hosts an intermediate-mass 
YSO driving a 4 M0 outfiow with a mechanical luminosity 
of 5.6 L0. The outfiow has been map ped with the IRAM 
30-m telescope in a number of species iCodella fc Bachillen 
I1999I) and four molecular clumps were identified along the 
main axis. The physica l conditions derived in the four 
clumps by ICodella fc Ba chiller ( 1999) by means of LVG or 
LTE calculations are very similar to each other, and we use 
these values to guide the choice of the input parameters in 
our chemical models. We choose to fit the clump SI of the 
CB3 outflow located in the southern lobe at an offset (0",- 
30") from the millimeter driving source since it is the best 
defined clump. In Table we list all the species and the 
lines observed in the SI position for which the molecular 
data (Einstein coefficients and coUisional rates) exist and 
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where i is the index indicating the line. A'' is the total num- 
ber of lines considered, FluXmod is the integrated ffux of the 
model and Fluxots is the observed integrated fiux. The dif- 
ference between the theoretical and observed fiux is weighted 
with the observed flux instead of the error of the measure- 
ment because the uncertainties associated to the chemical 
models are higher than the error of the observations, mak- 
ing it meaningless to use the error as the weight. Moreover, 
in this way we avoid giving the brighter lines a higher weight. 
A similar formula has also been used by other authors, e.g. 
|Dot£ei^(|200i). 



4 RESULTS 

In the next 5 subsections, we will analyse the results we 
obtain for the fit of each molecular species separately; in 
Section 5 we will attempt to fit all the species simultane- 
ously. In order to take into account the uncertainties in the 
parameters adopted in the models and the calibration errors 
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Table 2. List of the models and their parameters: presence (Y) or absence(N) of non-dissociative shock, density, final temperature, 
percentage of freeze-out, initial sulfur abundance and other possible parameters that are different with respect to the other models. For 
Grid D models the density and temperature are the final ones (at t=10"' yr). 



model 



Al 
A2 
A3 
A4 
A5 
A6 
A7 
A8 
A9 
AlO 



shock 



n(H2) 



N 
N 
N 
N 
N 
N 
N 
N 
N 
N 



(105 



T 
(K) 



FR 

(%) 



X(S) 
(10-^) 



1 
1 
1 
1 

10 
10 
10 
10 
50 
5 



210 
210 
210 
210 
210 
210 
210 
210 
210 
210 



15 
35 
11 
20 
20 
30 
55 
80 
25 
30 



130 
130 
130 
19 
130 
1.3 
1.3 
1.3 
130 
1.3 



Bl 


Y 


5-10 


110 


40 


1.3 


B2 


Y 


5-10 


110 


60 


1.3 


B3 


Y 


5-10 


110 


60 


1.3 


B4 


Y 


5-10 


110 


80 


1.3 


B5 


Y 


5-10 


110 


60 


1.3 


B6 


Y 


5-10 


110 


60 


1.3 


B7 


Y 


5-10 


210 


60 


1.3 


B8 


Y 


2.5-5 


210 


40 


1.3 


B9 


Y 


3-6 


210 


60 


1.3 


BIO 


Y 


5-10 


210 


20 


1.3 


Bll 


Y 


3-6 


210 


40 


1.3 


B12 


Y 


2.5-5 


210 


60 


1.3 


CI 


N 


5-10 


110 


40 


1.3 


C2 


N 


5-10 


110 


60 


1.3 


C3 


N 


5-10 


110 


60 


1.3 


C4 


N 


5-10 


210 


60 


1.3 


Dl 


N 


5-10 


110 


20 


1.3 


D3 


N 


5-10 


110 


20 


1.3 


D4 


N 


2.5-5 


100 


20 


1.3 


D4-shock 


Y 


5-10 


100 


20 


1.3 


D8 


N 


5-10 


100 


20 


1.3 



Note 



retarded collapse 
retarded collapse 



Asl 


Y 


5 


210 


30 


1.3 


As2 


Y 


1 


210 


30 


1.3 



100%H2CO^CH3OH 

lOOyoHaCO^CHgOH 

10%CO=>CH3OH 

20%CO^CH3OH 



100%H2CO=>CH3OH, time-dependent evaporation 



Voutflow= 2 km s 1 
Voutflow= 100 km s-l 

^outfloiu~ 
^outflow~ 



= 2 km s-l 

= 2 km s-l 

20 km s-l 



of the observations, we assume that we obtain a good match 
to the models when the theoretical line flux differs less than 
10% from the observed flux. 

4.1 SO 

As already noted in Paper I, if in the chemical model we as- 
sume a standard solar initial abundance of sulfur (1.3x10^'') 
the final abundance of SO is much higher than the observed 
one. Moreover the radiative transfer model produces strong 
self-absorption in the lines, which is not present in the ob- 
served line profile. The observed SO lines can be better re- 
produced assuming a depleted sulfur initial ab undance of 
1.3xl0~^ as also indica ted by previous studies l|Viti et alJ 
Ein^:lRnflReet~d1lT99l. 

We have detected two lines of SO and, due to the dif- 
ferent excitation temperature of these lines (24 and 11 K), 
their ratio should allow us to constrain the temperature of 
the gas. Indeed, we see that models with temperature of 
T~200 K are the only models that can reproduce the ob- 
served line ratio. In Tablej^lwe list all the models where both 
the two SO line fluxes differ ^10% from the observed fluxes. 



ordered by increasing x^- The SO column density in all the 
models ranges from 3.6xl0"cm"^ to 6.4xl0"cm"^. The 
models from Grid B are the best ones, while models with- 
out a non-dissociative shock and models with low density 
(n(Il2)<3xl05 cm"^) do not seem to flt the SO lines at aU. 
In Table^we note that models with equal input parameters 
but two different percentage of freeze-out (BIO and B7 or 
B9 and Bll) can both flt the observed lines but at different 
epoch. In particular, the models with the lower percentage of 
freeze-out fit the lines at an early epoch with respect to the 
models with higher freeze-out. The best fit is given by the 
model BIO (n(H2)=(5-10)xl05 cm-^ T~200 K, FR=20%, 
Shock=Y) at a time t=10^ yr with N(SO)~4xl0^*cm"^ 
However, in the BIO model the SO column density increases 
with time, it has a maximum at t=3xl0^ yr and then it 
decreases again to a value similar to the one assumed at 
t=10^ yr. This implies that a good fit can also be obtained 
at a later time of t=9 x 10^ yr. The comparison between the 
observed lines and the best fit model is shown in Fig. |21 
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Table 3. List of the best models for the two SO transitions, 
ordered by increasing x i with the relative epoch, gas density, 
temperature, freeze-out and SO column density. 



Table 4. List of the best models for the CS(3-2) transition, or- 
dered by increasing x^i with the relative epoch, gas density, tem- 
perature, freeze-out and CS column density. 



model 


epoch 


n{H2) 


T 


FR 


N(SO) 




(yr) 


(10^ cm-3) 


(K) 


% 


(cm-2) 


BIO 


1x103 


5-10 


210 


20 


3.8x101* 


AlO 


5x10^ 


5 


210 


30 


6.4x10^* 


Bll 


2x10^ 


3-6 


210 


40 


5.5x10" 


B9 


3x10* 


3-6 


210 


60 


5.3x10" 


BIO 


9x10* 


5-10 


210 


20 


3.8x10" 


B7 


8x10* 


5-10 


210 


60 


3.6x10" 



D,s - S0(5s-b, 



D.2 - 





Figure 2. The two SO 65 — 54 and 43 — 32 lines observed in the 
clump SI of the CB3 outflow (histogram) with the line profiles 
predicted by the best fit model BIO at t=lxlO'' yr (continuum 
line). 



4.2 CS 

CS is a high density tracer and we find that at a fixed den- 
sity it is also quite sensitive to the percentage of freeze-out, 
and in fact models with the same temperature and density 
but different freeze-out fit the line at different epochs. In 14 
out of 33 models, at least a time can be found in which the 
flux of the CS(3-2) line differs less than 10% from the ob- 
served flux (see Table |3|l and for 3 models, namely Asl, C3 
and B9, we find a range of times fitting the line. The mod- 
els cover almost all the explored parameter space indicating 
that the CS(3-2) line alone can not be used to constrain the 
physical conditions of the clump, nor can it be used to dis- 
criminate among different scenarios of formation. Despite 
the different physical parameters, all the models in Table 
0have a similar CS column density: N(CS)~10i*cm"^. It 
is worth noting that the best fit model is the model Asl 
(n=5xl0^ cm-^ T~200 K, FR=30%, Shock=Y) at a time 
t=2xl0^ yr with N(CS) = 1.2xlO"cm"^. The comparison 
between the observed line and the best fit model is shown 
in Fig. 13 The fact that the modelled line has a peak higher 
than the observed one suggests that either the line has a 
FWHM slightly lower than the fixed value of 10 km s"^ or 
the line is weakly self-absorbed. 



4.3 SO2 

In Table |S] we report all the models where the SO2 line flux 
differs SC10% from the observed flux, ordered, as usual, by 
increasing y^ . The SO2 column density ranges from 2x10^* 
to 5x10^* cm~^. We find that N(S02) is very sensitive to 
time and in all the models it increases by more than two 
orders of magnitude from t=10^ yr to t=10^ yr. However, 



model 


epoch 


n(H2) 


T 


FR 


N(OS) 




(yr) 


(lO^cm-^) 


(K) 


% 


(cm-2) 


Asl 


2x10^ 


5 


210 


30 


1.2x10" 


CI 


1x10^ 


5-10 


110 


40 


1.2x10" 


A4 


lxl03 


1 


170 


20 


1.5x10" 


A2 


3xl03 


1 


210 


35 


1.3x10" 


B12 


4x10* 


2.5-5 


210 


60 


1.4x10" 


BIO 


5xl03 


5-10 


210 


20 


1.6x10" 


D4 


2x10^ 


10 


70 


20 


2.2x10" 


Asl 


1x10^ 


5 


210 


30 


l.lxlO" 


03 


6x10^ 


5-10 


110 


60 


l.lxlO" 


D8 


2x10^ 


10 


70 


20 


1.6x10" 


B9 


4x10* 


3-6 


210 


60 


1.6x10" 


03 


5x10^ 


5-10 


110 


60 


1.2x10" 


B4 


3xl03 


5-10 


110 


80 


l.OxlO" 


02 


4x10* 


5-10 


110 


60 


l.OxlO" 


B9 


3x10* 


3-6 


210 


60 


1.3x10" 


Asl 


3x10^ 


5 


210 


30 


1.3x10" 


03 


5xl03 


5-10 


110 


60 


1.7x10" 


B8 


1x10* 


2.5-5 


210 


40 


l.OxlO" 


As2 


1x10^ 


1 


210 


30 


1.2x10" 



CS(3-2) 




Figure 3. The OS(3-2) line observed in the clump SI of the CB3 
outflow (histogram) and the line profile predicted by the best fit 
model Asl at t=2xl03 yr (continuum line). 



we do not seem to find a good flt for times > 10 yr. At later 
times the modeled column density is usually higher than the 
observed one. Similar to our SO fits, the models from grid B, 
As and A with high H2 density (nJSlO cm~ ) give the best 
fits, indicating once again that the scenario where the clump 
is (at least partially) pre-existing the outflow is the most 
likely. All models from Grid D, C and A (with the exception 
of A8 and A9 - see below) produce a SO2 line intensity 
very different from the observed one, up to two orders of 
magnitude for the model D4. One common characteristic 
of these models is the lack of a post-shock phase at high 
temperature, hence it seems that this species is indeed a 
good tracer of the presence of a non-dissociative shock. This 
is in agreement with theoreti cal models of C-type shocks 
JPineau des Forets et alJll993i) where the abundance of SO2 
is seen to increase after the passage of the shock and SO2 is 
the dominant form of sulfur in the post-shock region. 

The best fit model (see Fig. 2J is the model A9 
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Table 5. List of the best models for the SO2 3i3-2o2 transition, 
ordered by increasing x i with the relative epoch, gas density, 
temperature, freeze-out and SO2 column density. 



Table 6. List of the best models for the H2CO 321-220 transition, 
ordered by increasing x^i with the relative epoch, gas density, 
temperature, freeze-out and H2CO column density. 



model 


epoch 


n{H2) 


T 


FR 


N 


(SO2) 


model 


epoch 


n(H2) 


T 


FR 


N 


H2C0) 




(yr) 


(10^ cm-3) 


(K) 


% 


(cm-2) 




(y) 


(10^ cm--^) 


(K) 


% 


(cm-2) 


A9 


6x10^ 


50 


210 


25 


2.5x10^'' 


Asl 


1x10^ 


5 


210 


30 


6.1x10" 


Bl 


3x10=* 


5-10 


110 


40 


3.5x101'' 


B3 


1x10^ 


5-10 


110 


60 


6.5x10" 


A8 


2x10^ 


10 


210 


80 


2.6x101** 


Bll 


4x10^ 


3-6 


210 


40 


6.5x10" 


B8 


5x10^ 


2.5-5 


210 


40 


2.7x10" 


BIO 


4x10^ 


5-10 


210 


20 


5.1x10" 


A8 


3x10^ 


10 


210 


80 


2.7x10" 


Asl 


2xl03 


5 


210 


30 


6.6x10" 


B7 


6x10^ 


5-10 


210 


60 


5.1x10" 


B8 


3x10=* 


2.5-5 


210 


40 


6.0x10" 


Asl 
B9 


9x10^ 


5 


210 


30 


2.1x10" 


BIO 


3x10^ 


5-10 


210 


20 


4.5x10" 


6x10"^ 


3-6 


210 


60 


3.3x10*'* 














B3 


1x10=* 


5-10 


110 


60 


3.8x10" 














AS 


4x10^ 


10 


210 


80 


2.8x10" 














Bll 


7x10^ 


3-6 


210 


40 


3.5x10" 
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Figure 4. The SO2 3i3-2o2 line observed in the clump SI of the 
CB3 outflow (histogram) with the line profiles predicted by the 
best fit model A9 at t=6xl0^ yr (continuum line). 



Figure 5. The H2CO 32i-22o line observed in the clump SI of 
the CB3 outflow (histogram) with the line profiles predicted by 
the best fit model Asl at t=lxl0'^ yr (continuum line). 



(n=5xlO'' cm-^ T~200 K, FR=25%, Shock=N) at the 
epoch of t=6xl0^ yr. We do not think, however, that this 
model is correct because such high uniform density for a 
clump of 0.12 pc would imply a visual extinction of more 
than 1000 mags. Moreover, the initial sulfur abundance for 
A9 is solar. Probably, in the A9 model the absence of the 
non-dissociative shock is compensated by the higher initial 
sulfur abundance and the high density that lead to a SO2 
column density able to fit the line. As can be seen from Ta- 
ble 13 there are quite a number of models that can fit the 
SO2 line. It is clear therefore that we can not use this line 
to deduce uniquely the physical characteristics of our clump. 
However, the SO2 line does seem to trace a gas with a density 
between 3x10^ and 10® cm~^ and a temperature between 
100 and 200 K. The freeze-out parameter is not constrained 
at all since in the models that fit the data the freeze-out 
parameter ranges from 25% to 80%. This is not surprising 
since SO2 is a second generation species which forms after 
the dissociation of other sulfur bearing species released from 
the grains. 



4.4 H2CO 

In Table |S| we report all the models that fit the observed 
flux of the H2CO 321-220 line within 10%. Good fits can 



be obtained in models with non-dissociative shocks (Grid B 
and As). In the models without them (Grid A, C and D) 
the formaldehyde is usually overabundant. This is because 
it is mainly formed by the reaction of oxygen with CH3. In 
models with a non-dissociative shock, during the high tem- 
perature phase, the oxygen mainly interacts with neutrals 
to form OH and H2O, thus there is less oxygen free to form 
H2CO. H2CO is one of the species most sensitive to time; 
indeed its abundance increases by more than one order of 
magnitude from t=10'^ yr to t=10^ yr. A good fit can only 
be found for very early epochs of ^ 4 x 10^ yr. As in Pa- 
per I we find that if the clumps evolve for more that 5 x 10^ 
yr after the high temperature phase, the theoretical abun- 
dance of H2CO is always overabundant with respect to the 
observations. Even in models where we allow conversion of 
H2CO to methanol on the grains (B3 and B4) an acceptable 
fit can only be obtained at t=10'^ yr for the B3 model. The 
H2CO column density in the selected models ranges from 
4.5x10*'* to 6.5x10** cm~^ and the freeze-out parameter 
ranges between 20% and 60% while the temperature is ^^200 
K, with the exception of the B3 model. The best fit model 
is Asl (n(H2)=5xlO^ T~200 K, FR=30%, Shock = Y) at 
the time t=10^ yr and column density N(H2CO)=6.1xlO*'* 
cm~'^. In Fig. I^the Asl model is plotted superimposed to 
the observed line spectrum. 
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4.5 CH3OH 

Methanol has always been considered a Ro od tracer of the 
physical condition of the emitting gas (e.g. iKalenskii et alJ 
I1997I . iLeurini et al . 2004). In particular, it has bee n used as 
probe of density in high-density medium (e.g. Menten et alJ 
Il988r and a temperature estimate can be easily derived from 
the rotational diagram, if multiple transitions are available. 
Moreover, the presence of several transitions close in wave- 
length make it easy to observe more lines simultaneously, 
increasing the number of observational constraints and min- 
imizing the relative calibration uncertainties. 

The Boltzmann plot of the 8 methanol lines observed 
in the SI clump gives a column density N(CH30H)~ 
10^^ cm~^ and a rotatio nal temperature of only 16 K 
JCodella fc Bachilleilll999h . well below the 200 K we derived 
from the other molecules. However, in the physical condi- 
tions typical of outflows the gas is not in LTE condition and 
it is subthermally excited so that the rotational tempera- 
ture estimate d from the Boltzmann plot is below the kinetic 
temperature JBachiller et al.lll995l) . 

In the radiative transfer model the A and E types of 
methanol must be considered as separate species since the 
transition between the two symmetry states can not happen 
at a signiflcant rate by radiative or coUisional processes. We 
assume the standard ratio of 1 between A and E methanol 
as also suggested by the absence of a signiflcant shift in the 
A and E transitions in the Boltzmann plot. This choice has 
also been verified a posteriori since we find that both A and 
E methanol lines have similar deviation with respect to the 
observed lines. 

From Tabled we see that the 3k-2k and 2k-1k tran- 
sitions are very close in frequency. The line rest frequen- 
cies of the CH3OH 3o-2o A+ (145.103 GHz) and 3-i-2_i 
E (145.097 GHz) transitions are separated by 12.4 km s^^, 
those of the 3-i-2_i E (145.097 GHz) and 3o-2o E (145.094 
GHz) by 6.2 km s'^ . The 2o-loA+ (96.741 GHz) and 2_i- 
1-1 E (96.739 GHz) transitions are also separated by 6.2 
km s~^. Since the separation of the lines is similar or lower 
than the typical FWHM in the clump (10 km s~^) we expect 
that 3k-2k and 2k- Ik transitions overlap with each other. 
Indeed the lines at 145.1 and 96.7 GHz are blended; thus 
photons at these frequencies will see a higher optical depth 
than expected, because they are "seeing" the line profiles 
of neighbouring transitions. The SMMOL code does not yet 
have a line-overlap facility, so our modelling of these blended 
methanol transitions is semi-quantitative. For this reason we 
consider as an acceptable fit for the methanol all the mod- 
els that produce line fluxes that differ less than 30% with 
respect to the observed flux, instead of 10% as used for the 
other species. 

None of our models are able simultaneously to fit all of 
the 8 methanol transitions (5 lines of CH3OH-E and 3 lines 
of CH3OH-A). In particular, we find that in models where 
the higher energy transitions (5k-4k) are reproduced, the 
lower energy transitions (3k-2k and 2k- Ik) are largely un- 
derestimated. Conversely, models that reproduce the lower 
energy transitions overestimate the transitions at higher en- 
ergy. This could be a direct consequence of the effects of 
line-overlap in the 3k-2k and 2k- Ik transitions: the low en- 
ergy transitions are simply "seeing" more optical depth than 
the higher energy transitions. 



Table 7. List of the best models for the 3 Unes of the CH3OH 
5k-4k transition (both A and E type), ordered by increasing x'^, 
with the relative epoch, gas density, temperature, freeze-out and 
methanol column density. The last 2 models have a smaller clump 
size of 0.05 pc 



model epoch 

(yr) 



n(H2) T FR N(CH30H) 

(10^ cm-3) (K) % (cm-2) 



Dl 


(9-10) xlO"* 


5-10 


110 


20 


(4.2-4.8) xlO^"* 


B12 


(8-9)xl03 


2.5-5 


210 


60 


(4.3-4.7) X 10" 


D8 


1x10^ 


5-10 


100 


20 


5.2x10" 


D4 


7x10* 


3-5 


100 


20 


4.3x10" 


AlO 


(1-10) X 10-1 


5 


210 


30 


(1.2-1.5)xl0i5 


Asl 


(1-4) xlO^ 


5 


210 


30 


1.5xl0i5 



We did test, however, whether these inconsistencies are 
due to incorrectness of our basic assumptions on the phys- 
ical condition of the clump; we computed further models 
where we assumed: i) a lower temperature of T=20 K or ii) 
a smaller size of the clump of 0.06 pc. Again, there is no 
single model able to fit both the lower and the higher en- 
ergy transitions at the same time. Another possible expla- 
nation is that the high energy transitions trace a different 
component of the gas, possibly with different physical (and 
chemical) characteristics. Indeed the 5k-4k lines seem to be 
fitted with models that produce a methanol column den- 
sity N(CH3OH)~5xl0" cm"^ while the 3k-2k and 2k-1k 
lines are fitted by models that have a higher column density 
N(CH3OH)~5xl0^'^ cm"^. Hence, we tried to separately fit 
the two components. For the high energy transitions (the 
5k-4k lines) we find that most of the models that fit the 
lines belong to Grid D (see Table |7J, indicating that these 
transitions may trace a gas compressed by the outfiow (sce- 
nario D). Note that the models in Grid D have a lower abun- 
dance of methanol than the models in the other grids. This 
is due to reduced amount of frozen material, since methanol 
is preferentially formed on the grains via hydrogenation of a 
fraction of CO and H2CO: in the scenario depicted in Grid 
D the clump is formed by the compression of the outflow, 
and therefore the timescale available for freeze out is quite 
short, as explained in detail in Paper I. Moreover during the 
clump formation the temperature also increases, slowing sig- 
nificantly the freeze out. While it is possible that the high J 
transitions trace a newly formed clump (scenario of Grid D), 
it is also possible that this lower column density component 
simply traces a smaller component of the clump. To test this 
hypothesis we recomputed models AlO and Asl assuming a 
smaller clump size of 0.04, 0.05 and 0.06 pc instead of 0.12 
pc. We found that both models AlO and Asl were able to 
fit the 3 CH3OH lines of the 5k-4k transition within 30% 
with a clump size of 0.05 pc (see Table 0. 

As discussed above the three 3k-2k lines and the two 
2k- Ik lines are closely blended due to line-overlap effects. 
Therefore to attempt a fit to the 5 CH3OH lines of the 3k-2k 
and 2k- Ik transitions, we decided to calculate the x^ con- 
sidering the sum of the fiux of all the blended lines rather 
than each line fiux separately. We find that two of our mod- 
els fit the 5 lines (see Table |HJ. The best fit model is BIO 
(n(H2) = (5-10)xl0'5 cm-^ T~200 K, FR=20%, Shock = Y) 
with a column density N(CH3OH)=7xl0^^ cm"^ The best 
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Table 8. List of the best models for the 5 lines of the CHsOHSk- 
2k and 2k-1k transitions (both A and E type), ordered by in- 
creasing x^i with the relative epoch, gas density, temperature, 
freeze-out and methanol column density. 



model 


epoch 

(y) 


n{H2) 
(10^ cm-3) 


T 

(K) 


FR 

% 


N(CH30H) 
{cm-2) 


BIO 
C3 


6x10'' 
7x10* 


5-10 
5-10 


210 
110 


20 
60 


3x1015 




20- CHpH(2,-lJ 



0,0 ^^ ^u^^ 



i V^^ 



Figure 6. The CH3OH 3k-2k and 2k-1k lines observed in the 
clump SI of the CB3 outflow (histogram) and the line profile 
predicted by the best fit model BIO at t=6xl0* yr (continuum 
line) . 



fit time is t=6xl0* yr, although in BIO the methanol col- 
umn density is quite constant in time so it is not possible to 
constrain the age of this component. In Fig. El we show the 
best fit model: although the total flux of the three blended 
3k-2k hnes are still within the 30% of tolerance, we do not 
achieve a good fit for the CH3OH 3o-2o A"*" line at -41.5 km 
s~^. We note that in all the models the ratio between the 
3o-2o A^ and the 3o-2o E lines is always higher (up to a fac- 
tor of 1.5) than the observed value indicating a systematic 
trend to overestimate the 3o-2o A+ irrespective of the model 
parameters. This is due to the effect of line overlap that is 
not included in our radiative transfer model. 

Despite the limitation due to the lack of line-overlap 
effect in the modeling, our analysis seems to indicate the 
presence of two different methanol components: one with 
lower column density N(CH30H) ~ 5x10^-10^'^ cm'^ and 
possibly with smaller size, traced by the higher energy tran- 
sitions, and another with higher column density N(CH30H) 
~ 7x 10^5 cm^^ and probably at lower excitation conditions 
since it is traced by the lower energy lines. In our analysis 
we assumed the same temperature for both components; this 
may be unrealistic but we believe that in any case neither of 
the two components has a temperature less than 100 K, due 
to the presence of the outflow. The possibility of different 
components or substructures within the clump was already 
put forward in Paper I for CB3 and partly confirmed by the 
better fits obtained for L1157 where the substructure is in 
fact observationally resolved. 



5 DISCUSSION 

In the previous section we have shown that, if one con- 
siders each species separately, it is indeed possible to find 
a very good fit for each species. The column densities we 



find are always within a f actor of 3 of the ones derived by 
ICodella fc Bachilleil ^1999^ by using LTE or LVG models. 

We now attempt to find a model that simultaneously 
fits all the molecular transitions observed in the SI clump. 
In order to take into account the large uncertainties in the 
model for the global fit we consider acceptable all the models 
where the line flux differs less than 50% from the observed 
flux. We choose this relatively high level of tolerance be- 
cause the line emission depends, directly or indirectly, on 
many parameters that are not well known such as the size 
of the emitting region, the density, the temperature and the 
chemical abundances that in turn depend on other parame- 
ters such as the percentage of depletion of the gas on to the 
grains, the efficiency of some reactions, the density structure 
etc, and each species/transition has a different sensitivity to 
each parameter. We find that none of the 33 models is able 
to fit simultaneously all the molecular lines within 50% in- 
dicating, as also suggested in Paper I, that indeed several 
gas components exist in the IRAM-30m beam. 

In order to identify the different components we con- 
sider different groups of molecules. Firstly, we consider all 
the S-bearing molecules i.e. CS, SO, SO2. None of the 33 
models are able to fit all the 4 observed lines. In particu- 
lar, while several models fit both the CS and SO transitions 
(see Table |UJ the same models do not fit the SO2 line that 
is usually fitted at the earlier times of t<10* yr(see Table 
|3l. The best fit model for CS and SO is the model B9 at a 
time of 1=3x10* yr (see Fig.0; however, at that time, the 
S02(3i3-2o2) line is 8 times more intense than the observed 
emission (see Fig. |S| where the column density of the species 
is plotted versus time). In the same model the SO2 line can 
however be fitted at an early time, t=6xl0^ yr, when the 
column density is N(SO2)=3xl0i* cm~^. Since in most of 
the models the SO2 abundance increases significantly with 
time, the fact that the S02(3i3-2o2) line is fltted by young 
epoch (see Table I^J indicates that a low quantity of SO2 is 
needed to flt the line. A low column density of SO2 can be 
obtained not just with a low chemical abundance but also 
with a smaller size of the clump. In fact a smaller component 
in the clump has been already suggested by the analysis of 
the methanol lines. As we did for the methanol, we test this 
possibility considering the chemical models AlO and Asl 
with a smaller clump size of 0.04, 0.05 and 0.06 pc. In the 
AlO models the SO2 abundance is always lower than that in 
the Asl model due to the lack of a high temperature phase 
in the former model. The AlO model does not succeed in 
fitting the observed S02(3i3-2o2) line. On the contrary, a 
very good fit can be obtained if we assume a size of 0.05 
pc for the model Asl at the epoch 1=4x10* yr. At this 
epoch, the model Asl also fits the 3 higher energy lines of 
methanol (see Fig.|^, as already shown in Section 4.5. The 
same model does not fit the CS and SO lines. We therefore 
conclude that SO2 and the higher energy methanol 5k-4k 
lines may be emitted from a smaller gas component with a 
size of d~0.05 pc, n(H2)- 5 x 10^ cm"^ and T~200 K. This 
component is associated with emission at high excitation of 
molecules produced in high temperature gas thus it could 
trace the zone of stronger interaction between the outfiow 
and the clump. 

Next we attempt to fit simultaneously the CS and SO 
transitions together with the 5 lower energy transitions of 
CH3OH (3k-2k and 2k- Ik). We find no models able to fit 
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Table 9. List of the best models for the SO(65 - 
and CS(3-2) lines, ordered by increasing x^i 



■54), SO{43-32) 
with the relative 
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Figure 7. The CS and SO lines observed in the clump SI of the 
CBS outflow (histogram) with the line profile predicted by the 
best fit model B9 at t=3xl0* yr (continuum line). 




Figure 8. Column densities versus time for the B9 model. 
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Figure 9. The S02(3i3-2o2) and CH3OH (5k-4k) lines observed 
in the clump SI of the CB3 outflow (histogram) with the line pro- 
file predicted by the best fit model Asl at t=4xl0^ yr (continuum 
line). The size of the clump is d=0.05pc. 



epoch, X I gas density temperature and freeze-out. 
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n(H2) 


T 
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(10^ cm-3) 


(K) 
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B9 


3x10^ 


3-6 


210 


60 


Asl 


1x10^ 


5 


210 


30 


B4 


2xl03 


5-10 


110 


80 


B12 


2x10" 


2.5-5 


210 


60 


B12 


SxlO'' 


2.5-5 


210 


60 


Asl 


2xl0» 


5 


210 


30 


Al 


2x10^ 


1 


210 
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simultaneously the 3 species within 50%. In fact, from the 
fitting of the single species (see Section 4.5) we see that the 
low energy methanol lines seem to trace a slightly higher 
density (model BIO with n(H2)=(5-10)xlO'^ cm"^) than the 
CS and SO component (model B9 with n(H2)=(3-6)xlO^ 
cm~ ) , but the same temperature (T~200K) and very simi- 
lar time (6x10'* yr and 3x10'' yr respectively). However the 
difference in the density is not so relevant, considering also 
the high uncertainties in the model especially for methanol 
because of line overlap effects, and we think that probably 
the low energy transitions of CH3OH could be associated 
with the same extended gas component traced by CS and 
SO. 

Finally we check whether the H2CO 32i-22o line could 
be associated with one of the two identified components. 
We find that the B9 model with clump size of 0.12 pc at 
t=3xl0* yr predicts a line about 25 times brighter than 
the observed one. Indeed, as said in Section 4.4, cissuming 
a size of 0.12 pc the H2CO can be fitted only at very early 
epoch (t^4xl0^ yr) while for more evolved times all the 
models predict far too much II2CO. On the other hand, if 
we consider the smaller component, i.e. the Asl model at 
4x 10"* yr with clump size of 0.05 pc, the theoretical intensity 
is only 3 times higher than the observed one, while at the 
early epoch of 2 x 10* a very good fit (within 20%) of the 
observed line is obtained. Although a global fit of the H2CO 
with the SO2 and high J CH3OII lines cannot be obtained 
within the 50% of tolerance, we think that H2CO does not 
trace the whole clump but is associated with the smaller 
gas component and it is probably emitted by a region even 
smaller than 0.05 pc. 

The SI clump has also been observed in the CO (2-1) 
line; however the profile of this transition is very different 
with respect to all the other lines: a blue non-Gaussian line 
wing and a self absorption at the ambient LSR velocity (- 
38.5 km"'^) are present in the line. This profile indicates that 
the line has a substantial contribution from the large scale 
outfiowing gas and also that part of th e emission is self ab- 
sorbed by the ambient cloud. Moreover. ICodella fc BachilleJ 
ilQQgtl estimate a temperature of ~10 K for the CO at all the 
line velocities and they find that the (2-1) line is optically 
thick showing that the emission is dominated by the outer, 
colder part of the outflow rather than by the SI clump. All 
these considerations lead us to not consider the CO (2-1) 
line for the global fit of the SI clump. However, we still ran 
the radiative transfer code for CO in order to roughly evalu- 
ate if in the best fit models the contribution of the clump to 
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this line is negligible as we expected. Indeed, we find that in 
the B9 model, that fits CS and SO, the predicted CO(2-l) 
flux is only 23% of the observed flux. 

We find that some species are very sensitive to 
time, such as S-beari ng molecules, in agreement with sev- 
eral previous studies llHatchell et aLlllQQSLlViti et alJl200ll 
IWakelam et al.l 12001 ICodella et al.l l2005ll . and they can 
therefore be used to determine the age of the SI clump. 
We fit both the components of the clump with the same 
time ~ 3-4x10'' yr which is similar to the dynamical time 
of the SI c lump estimated from the SiO emission (t^ 10^ yr, 
ICodella fc Bachillcr 1999). Note however that our starting 
point (t=0) in all scenarios is set to be when the outfiow 
starts interacting with the surrounding environment hence 
the age we derive is not necessarily the age of the outfiow. 



6 CONCLUSIONS 

A time-dependent chemical model has been coupled with a 
radiative transfer code to model the observed molecular lines 
emission in the small chemically rich clump SI along the 
south lobe of the CBS outfiow. This proved to be a powerful 
tool in constraining the physical and chemical parameters of 
the gas, especially when there is more than one transition 
of the same molecule. 

We find that more than one model is able to fit the line 
emission of single species. Hence to constrain the physical 
parameters of our clump a multi-line analysis is needed. 

Different scenarios for the formation of the clump have 
been investigated within a large range of physical conditions. 
Our results show that most of the observed lines can not be 
fitted by models from the scenarios A, C and D. On the 
contrary, the models from scenario B or As seem to repro- 
duce most of our data. This is consistent with our previous 
findings (see Paper I) . We conclude that the observed dense 
clumps are, at least partially, formed prior to the advent of 
the outflow. The advent of the outfiow on such clumps leads 
to a short phase of high temperature (T^IOOOK) followed 
by fast cooling down to temperatures of ^ 200 K. 

In addition to the large scale outflowing gas traced by 
the non-Gaussian wings detected in the CO (2-1) line, we 
identify two other components inside the large IRAM-30m 
beam. The extended component (d=0.12 pc) is traced by CS 
and SO and probably by the lower energy transitions (3k-2k 
and 2k-1k) of CH3OH. This component is not uniform in 
density but it has a density structure in the range 3x10^- 
10^ cm~^. The second component is more compact, with 
size ~0.05 pc and density ~ 5x10^ cm~"^, and it is traced 
by SO2, H2CO and the 5k-4k CH3OH lines. For both the 
components we estimated a similar temperature (T~200 K) 
and time (t^SxlO'' yr) from when the clump has started to 
interact with the outflow. Since the compact component is 
traced by the higher energy transitions of CH3OH, it is prob- 
ably at higher excitation and may trace the zone of stronger 
interaction between the outflow and the clump. However, 
the analysis of the CH3OH data needs to be confirmed by 
modelling of the data with a radiative transfer code that 
takes into account the effects of line-overlap (Gray & Yates, 
in prep). 
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